Method of layer stripping to predict subsurface stress regimes

ABSTRACT

A method for analyzing seismic shear wave data, using a layer stripping technique, to predict subsurface stress regimes is disclosed. Polarization directions of shear wave data, from either a vertical seismic profile or from surface reflection data, are analyzed, and time lags between fast and slow split shear wave are determined. Natural polarization directions of and time lags between the split shear waves in an upper layer are determined above the shallowest depth where data cues suggest polarization changes take place. Source and receiver axes of the data below the depth of polarization changes are rotated by an azimuth angle, to bring the axes into proper aligment. A static time shift is then applied to eliminate the time lag in the upper layer above the depth where polarization changes were indicated.

FIELD OF THE INVENTION

The present invention relates generally to geophysical production of oiland gas. More specifically, this invention provides a method forreliably and accurately applying a layer stripping technique to predictsubsurface stress regimes.

BACKGROUND OF THE INVENTION

Shear wave (S-Wave) seismic exploration techniques have historicallyemployed shear wave seismic sources and shear wave seismic receivers ina seismic survey to gather seismic data. Such a seismic survey has beeneither linear or areal in its extent. The seismic energy imparted by theshear wave seismic source is detected by the shear wave seismicreceivers after interacting with the earth's subterranean formations.Such seismic surveys, however, until recently have been limited toutilizing a shear wave seismic source having a single line of action orpolarization, oriented with respect to the seismic survey line ofprofile, to preferentially generate seismic waves of known orientation,e.g., horizontal shear (SH) waves or vertical shear (SV) waves. Theshear wave seismic receivers utilized in conjunction with a given shearwave seismic source have similarly been limited to a single line ofaction or polarization, oriented with respect to the seismic survey lineof profile, to preferentially receive a single component of the seismicwave, e.g., (SH) wave or (SV) wave. As used herein, the term "line ofaction" generally comprehends a defined vector displacement, such as theparticle motion of the seismic wave. In present shear wave seismicsurveys, the lines of action of the seismic source and the seismicreceivers usually have the same orientation relative to the line ofprofile and if so are said to be "matched".

The term "polarization" in the context of seismic waves refers to theshape and spatial orientation of particle trajectories. Here we restrictthe term to mean only the spatial orientation of the line along which aparticle moves in a linearly polarized wave. Hence "polarization" and"polarization direction", as used here, both imply the spatialorientation of such a line, the latter term emphasizing the restrictionto linear rather than more general (e.g., elliptical) motion. A"polarization change", then, does not mean a change, for example, fromlinear to elliptical motion nor a polarity reversal but only a change inthe spatial orientation of the line along which a particle moves.

As long as seismic surveys were limited to seismic sources and seismicreceivers having a compressional (P) wave lines of action, satisfactoryresults were generally obtained irrespective of the orientation of theseismic survey line of profile with respect to the underlying geologicalcharacter of the subterranean formations. However, when the seismicsources and seismic receivers are of the shear wave type, i.e., eitherhorizontal shear (SH) wave or vertical shear (SV) wave, the orientationof the seismic survey line of profile and/or the line of action of theshear wave seismic source with respect to the geological character ofthe subterranean formations can determine whether or not meaningfulseismic data is obtained.

As understood by those skilled in the art, compressional (P) waves arelongitudinal waves where the particle motion is in the direction ofpropagation. Shear waves are transverse waves where the particle motionis in a transverse plane perpendicular to the direction of propagation.Two special classes of shear waves are defined herein. Specifically,horizontal shear (SH) waves where the particle motion in the transverseplane is further restricted to be perpendicular to the line of profileof the seismic survey (i.e., horizontal) and vertical shear (SV) waveswhere the particle motion in the transverse plane is further restrictedto be perpendicular to the horizontal shear (SH) particle motion.

As the orientation of the seismic survey line of profile is dependent onthe geological character of the subterranean formation, when matchedshear wave seismic sources and shear wave seismic receivers are used, itis known by those skilled in the art that shear wave seismic surveys areadversely affected by azimuthally anisotropic subterranean formations.Azimuthally anisotropic subterranean formations are likely to havevertical planes of symmetry. Because shear wave behavior is complicatedand generally uninterpretable when the symmetry planes are neitherparallel to nor perpendicular to the line of action of the shear wave,care must be taken to ensure that the seismic survey line of profile islaid out either parallel or perpendicular to the symmetry planes.

When the seismic survey line of profile is laid out either parallel orperpendicular to the symmetry planes, the utilization of matched sets of(SH) wave and (SV) wave seismic receivers and seismic sources haveprovided useful information regarding the geological character of asubterranean formation. Such a technique requires prior knowledge of theseismic velocity anisotropy of the subterranean formation to besuccessful.

The interaction differences of (SH) waves and (SV) waves have beenutilized to detect and measure the anisotropic properties of anazimuthally anisotropic subterranean formation when the seismic lines ofprofile are properly oriented with respect to the symmetry planes andmatched sets of shear wave seismic sources and shear wave seismicreceivers have been deployed in the seismic survey. In suchapplications, (SH) and (SV) shear wave seismic sources and seismicreceivers are utilized, but only in matched sets, i.e., (SH) shear waveseismic sources with (SH) shear wave seismic receivers and (SV) shearwave seismic sources with (SV) shear wave seismic receivers. However, ifthe seismic survey line of profile is not properly oriented with respectto the planes of symmetry, the seismic information observed can bedifficult to interpret at best.

The orientation of the seismic survey line of profile with respect tothe symmetry planes is critical. Consequently, utilization of matchedsets of shear wave seismic sources and shear wave seismic receivers haveproduced inconsistent results when the seismic survey line of profilehas not been properly laid out with respect to the anisotropicgeological character of the subterranean formations.

Those acquainted with the art of seismic exploration, especially inseismically virgin territory, realized that prior knowledge of thegeological character of the subterranean formations is generally notavailable prior to seismic exploration. The method and system ofgeophysical exploration of the present invention can be advantageouslyemployed without regard to or knowledge of the geological character ofthe subterranean formations and still obtain meaningful seismic data.

U.S. Pat. No. 3,302,164 relates to seismic exploration for detectingfluids in formations by obtaining a ratio of the velocities of shearwaves and compressional waves along a seismic line of profile. In orderfor the ratio to be obtained, however, the frequency spectra of thewaves introduced by a seismic source had to be controlled according tothe average velocity ratio expected to be encountered. An article,"Combined Use of Reflected P and SH Waves in Geothermal ReservoirExploration," Transactions of Geothermal Resources Council, Volume 1,May 1977, discussed tests made using both compressional and shear wavesin exploring for and evaluating geothermal reservoirs. U.S Pat. No.4,286,332 relates to a technique of propagating seismic shear waves intothe earth from compressional wave producing vibrators. U.S. Pat. No.4,242,742 describes a technique of obtaining shear wave seismic datafrom surveys where impact devices for waves are used as a seismic energysource.

S-wave birefringence, a property of elastic waves in anisotropic solids,is common for S-waves traveling vertically in crustal rocks. Earlymodels of anisotropic sedimentary rocks proposed by explorationgeophysicists were often transversely isotropic with verticalinfinite-fold symmetry axes. Such solids are not birefringent forS-waves with vertical raypaths Earthquake seismologists (e.g., Ando etal., 1983; Booth et al., 1985), however, found near-vertical S-wavebirefringence in earthquake data in the early 1980s. At the same time,oil companies recording three-component (3-C) seismic data independentlyfound vertical birefringence in hydrocarbon-bearing sedimentary basins.(Winterstein). Researchers from Amoco, Exxon, Chevron and ColoradoSchool of Mines documented this vertical birefringence for the firsttime publicly in 1986 at annual meetings of the EAEG and SEG (e.g.,Alford, 1986; Willis et al., 1986; Becker and Perelberg, 1986; Frasierand Winterstein, 1986; Martin et al., 1986). Since then much additionalevidence for vertical birefringence in sedimentary basins hasaccumulated (e.g., Squires et al., 1989).

A common model for vertical S-wave birefringence is extensive dilatancyanisotropy (EDA) proposed by Crampin et al (1984). The essential featureof this model is that horizontal stresses such as those from platetectonics create vertically oriented, fluid filled cracks or microcrackswhich cause anisotropy that, unlike transverse isotropy with a verticalaxis, will cause vertical S-wave birefringence. The validity of EDA asan explanation for vertical birefringence is not established, but it andvariants of it have proved useful as a framework within which to recordand interpret experimental data. An alternate model, which we call theNur model (Nur, 1971; Nur and Simmons, 1969), proposes the unstressedrock is isotropic with a uniform distribution of randomly orientedcracks. Axial stresses preferentially close the cracks perpendicular tostress directions, making the rock anisotropic. It is almost certain,whatever the best model proves to be, that much of the observed verticalS-wave birefringence results in some way from horizontal stresses.Crampin and Bush (1986) also pointed out that vertical S-wavebirefringence might provide a useful tool for reservoir development. Thepolarization direction of the fast S-wave in simple cases gives thedirection of maximum horizontal compressive stress, a quantity much indemand by those who induce fractures in reservoirs by techniques such ashydraulic fracturing.

Available evidence, (discussed later), including offset VSP informationsupports the notion that the vertical S-wave birefringence is caused byhorizontal stresses, and that the polarization direction of the fastS-wave lies in the direction of maximum horizontal compressive stress,even when subsurface structures are steeply dipping. It is likelyhowever that rocks exist for which the polarization direction of thefast S-wave for vertical travel does not lie along the maximumhorizontal stress direction. Rocks with fractures oriented by ancientstress regimes, or rocks of low symmetry with tilted symmetry axes, forexample, might constrain the fast S-wave polarization to lie in adirection other than that of maximum horizontal stress.

Unmistakable evidence is hereby presented for major changes in S-wavepolarization direction with depth (see also Lee, 1988). A relationshipbetween these polarization changes and any change of horizontal stressdirection certainly exists, and the S-wave birefringence data providepotentially useful information for reservoir development regardless whatthe relationship is. U.S. Pat. Nos. 4,803,666 and 4,817,061 (both toAlford) are hereby incorporated by reference. Alford discloses a methodof determining the S-wave polarization angles by finding the angle atwhich S-wave energy on off-diagonal components of an S-wave data matrixwas at a minimum. One implementation of Alford's method involvesselecting time windows that include only the leading portions of thefirst arrival S-waves, and then calculating energy on the off-diagonalcomponents at rotation angle increments of one degree.

However, an invalid assumption of Alford's rotation method is thatS-wave polarizations along a given raypath are generally orthogonal.Such an assumption is strictly valid only in certain symmetrydirections. The effectiveness of Alford's method is hindered by noise orby distortion of the signal on the off-diagonal components of the S-wavedata matrix.

Accuracy of analysis by Alford's rotation method depends, at least inprinciple, on having signal amplitudes of off-diagonal XY and YXcomponents identical at common times. If they are not identical, thedata do not fit the model, and the matrix cannot be diagonalized by asingle rotation of source and receiver coordinate frames. If signal onXY components differs systematically from that on YX components, therewill be systematic errors in calculated azimuth angles. But changes ofpolarization with depth cause just such systematic differences in signalon XY and YX components; specifically, the signal on one of the twocomponents lags that on the other by the amount imposed by the upperlayer.

Lefeuvre et al. (1989) and Cox et al. (1989) used propagator matrices ortransfer functions to analyze variations in S-wave birefringence withdepth in multicomponent VSP data, instead of applicant's proposed methodof layer stripping. These prior works utilize only a Fourier spectrum asan analytical method. Therefore, improvements in the S-wave data cannotbe readily seen, and the quality of the improvements do not matchapplicant's results. Being able to see the improved wavelet (as withapplicant's method) provides confidence to the analyst, as it providesinformation on how well the process is working.

Martin et al. (1986) analyzed changes in S-wave birefringence with depthin S-wave surface reflection data via a rudimentary layer strippingtechnique. They subtracted the effects of an upper layer to see theresidual effects in a lower layer. Their approach, however, required thegenerally unwarranted assumption that symmetry planes in a deeper layerwere orthogonal to those in an upper layer. That is, they did notperform any analysis to determine the actual orientation of the deepersymmetry planes.

Current methods of predicting subsurface fracture orientation or stressregimes fall short of providing accurate results, for the many reasonsdescribed above. There is therefore a need for an improved seismicmethod to evaluate changes in shear wave polarization with depth.

SUMMARY OF THE INVENTION

The present invention has been surprisingly successful in improving theanalyses of seismic shear wave data to predict subsurface stressregimes. Vertical seismic profile shear wave data or surface seismicreflection shear wave data has at least two linearly independent, nearlyorthogonal, and nearly horizontal source axes. Each source axis has atleast two corresponding receiver axes. An initial analysis of shear wavepolarization directions relative to a fixed coordinate frame is thenperformed, and apparent time lags between fast and slow shear waves aredetermined at several depths. Cues in the data are identified thatsuggest shear wave polarization changes.

The natural polarization directions of and the time lag between the fastand slow shear waves in an upper layer are determined, above andadjacent to the shallowest depth where the cues suggest polarizationchanges. Other depths may be used as well, even if no cues suggestpolarization changes The source and receiver axes of all the data thatare below or at the shallowest depth of indicated polarization changesare then rotated by an azimuth angle determined down to this depth, sothat the first source and receiver axes are aligned with the naturalpolarization direction of the fast shear wave, and the second receiveraxis is at a significantly different azimuth angle, and so that if thereis a second source, the second source and first corresponding receiveraxis are aligned with the natural polarization direction of the slowshear wave in the upper layer, while the second corresponding receiveraxis is at a significantly different azimuth angle.

A static shift is then applied to all data components corresponding toone of the effective sources, either to components corresponding to thesource aligned with the fast shear wave polarization direction, or tocomponents corresponding to the source aligned with the slow shear wavepolarization direction, to eliminate the time lag in the upper layerabove and adjacent to the shallowest depth where the cues suggestpolarization changes are indicated or suspected.

The invention may also be used for vertical seismic profile (VSP) dataor surface seismic reflection data that has only a single source axis.Only the receiver axes are rotated in this case.

If surface seismic reflection shear wave data is analyzed, one variationof the disclosed method includes an initial analysis of shear wavepolarization directions relative to a fixed coordinate frame insimilarly recorded VSP data from a nearby well, and the subsequentdetermination of the time lags.

A further variation of the invention permits analysis of surface seismicreflection shear wave data without the use of VSP data.

DETAILED DESCRIPTION OF THE DRAWINGS

FIG. 1 is a sectional view of the earth, illustrating the basic modelfor VSP shear wave recording.

FIG. 1a is a sectional view of the earth, illustrating the naturalcoordinate frame for vertical shear waves.

FIG. 2 is a sectional view of the earth illustrating the basic layerstripping rationale.

FIG. 3 is a plan view of the earth, illustrating the coordinate framefor recording and processing shear wave data, and the meaning of the 2×2shear wave matrix.

FIG. 4 shows the four shear wave components from the 1720 ft level ofwell 11-10X.

FIG. 5 shows the four shear wave components of FIG. 4 after "rotation".

FIG. 6 shows shear wave data from well 1-9J after "rotation".

FIG. 7 is a chart that illustrates polarization azimuths of the fastshear waves before layer stripping at the 1-9J well.

FIG. 8 is a chart that illustrates the polarization azimuths of the fastshear wave of the 1-9J well after layer stripping.

FIG. 9 is a chart that illustrates the polarization azimuths of the fastshear wave of the 1-9J well as a function of the initial rotation angle.

FIG. 10 is a chart that illustrates variations in shear wave lags withdepth, at the 1-9J well, after stripping off the near surface layer.

FIG. 11 is a chart that shows a summary of polarization angles of thefast shear waves with depth, for two independent layer strippinganalyses of the 1-9J well VSP data.

FIG. 12 is a chart that shows shear wave lag with depth, for the layerstripping sequence indicated by circles in FIG. 11.

FIG. 13 compares off-diagonal components of the 2×2 shear wave datamaterial of the 1-9J well before and after layer stripping.

DETAILED DESCRIPTION OF THE INVENTION

In accordance with the present invention, a new improved method andmeans for using layer stripping to predict subsurface stress regimes hasbeen developed.

The objective of the data analysis described herein is to quantifysubsurface shear wave (or S-wave) birefringence or, in other words, tofind the natural polarization directions of the two S-waves and the timedelays or lags between them. Natural polarization directions aredirections along which anisotropic rocks constrain polarizations ofS-waves to lie. The purpose of the analysis is to correlatebirefringence effects with formation properties such as direction ofmaximum horizontal stress. FIGS. 1 and 1a illustrate the basic model insimplest terms. An arbitrarily oriented horizontal displacement from asurface source propagates in the vertical direction as a fast S-wave(S₁), and a slow S-wave (S₂), with S₁ polarized along the direction ofmaximum horizontal compressive stress.

The term "polarization" in the context of seismic waves refers to theshape and spatial orientation of particle trajectories. The term isrestricted to mean only the spatial orientation of the line along whicha particle moves in a linearly polarized wave. Hence "polarization" and"polarization direction", as used here, both imply the spatialorientation of such a line, the latter term emphasizing the restrictionto linear rather than more general (e.g. elliptical) motion. A"polarization change", then, does not mean a change, for example, fromlinear to elliptical motion nor a polarity reversal but only a change inthe spatial orientation of the line along which a particle moves.

For arbitrary ray directions in anisotropic rocks of low symmetry, agreat deal of information is needed to interpret S-wave time lags andpolarizations. However, if the rocks have vertical twofold symmetryaxes, analysis is straightforward if raypaths are vertical, andpolarization directions relate in simple ways to symmetries of therocks. An initial assumption is that the rocks have vertical symmetryaxes and that their symmetry properties do not change with depth. Hence,in order to have raypaths as close to the symmetry axis as possible, thenear offset sources are positioned as close to the wells as possibleConcentric rings of offset VSPs serve primarily as a check on ourassumption of a vertical symmetry axis. That is, modeling showed that ifthe vertical direction is not a symmetry axis, S-wave polarizations atsmall offsets can vary asymmetrically with azimuth if the rocks are oforthorhombic or lower symmetry, even if there is a set of verticalcracks. On the other hand, if there is a vertical twofold symmetry axis,such S-wave polarizations will have twofold symmetry.

To determine natural polarization directions of the subsurface rock,several different rotation methods can be applied, as well as hodogramanalyses. The most reliable method in our experience is to find theangle at which S-wave energy on off-diagonal components of the 2×2S-wave data matrix is a minimum, a method we call the "rotation" methoddeveloped by Alford. All other methods had significant deficiencies. The"rotation" method can be implemented by choosing time windows thatinclude only the leading portion of the first arrival S-waves and thencalculating energy (sums of squares of amplitudes) on the off-diagonalcomponents at rotation angle increments of one degree. Only the leadingportions of wavelets need to be included because earlier observationsshowed that, after rotation to the angle which minimized off-diagonalenergy, the codas of diagonal wavelets differed from one another muchmore than did their leading edges. Hence, the leading edges are muchmore interpretable than the codas. The use of time windows provides aconsiderable signal-to-noise ratio (S/N) advantage over methods whichcalculate from individual points, and lends stability and consistency tothe answers. In most cases results are insensitive, within limits, tothe length of the time window.

An assumption of the "rotation" method, generally not valid, is thatS-wave polarizations are orthogonal. However, the assumption is strictlyvalid along any twofold symmetry axis and is a good approximation closeto such an axis. Seismic sources can be rotated by the same angle asreceivers, which is appropriate for vertical rays along verticalsymmetry axes in homogeneous anisotropic media. The differences inarrival times of fast and slow S-waves (the lags) can be computed bycrosscorrelating waves on the 2×2 S-wave matrix diagonals after rotatingto the angle that minimized off-diagonal energy. Lag is observed toincrease linearly with depth in a homogeneous, birefringent rock.

S-wave polarization directions were expected to remain constant withdepth, but data analysis showed convincingly that they did not.Polarizations at Lost Hills field changed relatively little; and if wehad considered only Lost Hills data, we probably would not have deemedit necessary to deal with polarization changes with depth. Polarizationchanges in Cymric and Railroad Gap fields to the south, however, werelarge and unmistakable, and a layer stripping method developed for datafrom those areas proved useful also for Lost Hills data.

Layer stripping involves simply subtracting off anisotropy effects in alayer in order to analyze anisotropy effects in the layer immediatelybelow. That is, S-wave splitting is cumulative, so that if anisotropychanges with depth, effects of anisotropy above the change, unlessremoved, will persist in the changed region and will confuse an analysisthere. Although polarization will change instantly when a wave enters aregion with different natural polarization directions, recorded waveletshapes change slowly and preserve information about their past travelsthrough other regions. Hence, if in polarization analysis one uses asignificant fraction of an arriving wavelet, as is done here, ratherthan just its "first arrival", which no one can accurately pick in realdata, one sees the effects of present as well as past polarizations.

What specifically hurts the effectiveness of the "rotation" 11 methodbelow a polarization change is distortion of signal on the off-diagonalcomponents of the 2×2 S-wave data matrix. Accuracy of analysis by the"rotation" method depends, at least in principle, on having signalamplitudes of off-diagonal XY and YX components identical at commontimes. If they are not identical, the data do not fit the model, and thematrix cannot be diagonalized by a single rotation of source andreceiver coordinate frames. If the signal on XY components differssystematically from that on YX components, there will be systematicerrors in calculated azimuth angles. But changes of polarization withdepth cause just such systematic differences in signal on XY and YXcomponents; specifically, the signal on one of the two components lagsthat on the other component by the amount imposed by the upper layer.This point can be understood by visualizing how wavelets in the lowerlayer of FIG. 2 project onto natural coordinate axes of the upper layer.

The inventive layer stripping process assumes certain subsurfaceproperties. For example, S-wave polarizations must remain practicallyconstant in a given layer. Polarizations hence are assumed to changediscontinuously at layer boundaries, and time lag in a given layerincreases monotonically from zero at the upper boundary to some finitevalue at the lower boundary. If polarizations were to changecontinuously with depth, the meaning of polarization analyses afterlayer stripping would be unclear. Also, each layer must be thick enough,and its birefringence large enough, to determine the correctpolarization direction and maximum lag for that layer. In ourimplementation, wave propagation is assumed to be close enough to asymmetry direction in every layer so that rotation of sources andreceivers by a single angle can do a good job of diagonalizing the 2×2S-wave data matrix.

To subtract off effects from above the depth at which polarizationchange occurs, all the data from below that depth is rotated by theazimuth angle determined down to that depth and then a static shift isapplied to remove the time lag between the two S-waves at that depth, asshown in FIG. 2. The S₁ and S₂ waves of the upper layer of FIG. 2 willact as independent sources, generating two sets of S₁ ' and S₂ ' wavesat the interface. Layer stripping removes the time delay between the twoeffective sources at the interface, causing the primed waves (S₁ ' andS₂ ') to behave as if the interface had been at the surface. The processsimulates putting a source at the depth where the polarization changeoccurs, such that the simulated source polarizations are oriented alongnatural polarization directions (assumed orthogonal) of the uppermedium. After layer stripping, rotation analysis is repeated as before,and further layer stripping (i.e., "downward continuation") is appliedif, for example, cues in the data indicate further polarization changes.

These layer stripping principles apply equally to surface seismicreflection data, but layer stripping will be less effective withreflection data because (1) signal-to-noise ratios are lower than indirect arrival VSP data, and (2) reflection events, which the methodmust rely on, do not necessarily occur close to where polarizationchanges occur. It may often be necessary to use information from VSPs tolayer strip surface seismic data.

Layer stripping, in contrast to methods involving the calculation ofpropagator matrices or transfer functions from depth to depth, typicallyexpects the user to judge where to do the stripping on the basis of apreconceived model; that is, he should have criteria in mind for judgingfrom analysis results where polarization directions change. Despite themore subjective nature of layer stripping vis a vis calculating transferfunctions, several possible advantages exist. First, layer strippingkeeps the user's focus on the geophysical objectives rather than detailsof calculations. Second, the user is able visually to evaluate effectsof stripping over large blocks of levels; this enables him to identifytrends and changes in trends without extra effort and thereby to picklayer boundaries perspicaciously. Third, layer stripping improves thequality of data for general interpretation.

It is usually necessary in any case to treat data in blocks of severallevels at a time, because it is impossible to determine birefringenceeffects if the two S-waves have not traveled long enough in thebirefringent medium to have accumulated a significant difference inarrival times. In noisy data, the robustness of birefringence analysisis aided by large lags between S-waves.

Cues that S-wave polarization directions have changed manifestthemselves as persistent changes with depth, in either the azimuthangles or the rate of change in time lags. Calculated azimuth anglestend to be insensitive indicators of polarization change below a thick,birefringent layer because properties of the S-wave wavelets remain muchthe same below the interface as they were above it, and the angles fromrotations consequently tend to remain the same for some distance belowthe change. In other words, S-wave splitting generates a kind of inertiain azimuth angle determinations. Lags, in contrast, are often sensitiveindicators of change: If polarization direction changes, the rate ofincrease in lag usually changes abruptly, and thus serves as theinterpreter's principal indicator of polarization change.

The procedure for layer stripping under normal circumstances may bedescribed in the following manner. The first step is to rotate sourceand receiver axes, say the x-axes, into alignment with the naturalpolarization direction of the fast S-wave in the upper layer. Therotation is applied to all data at and below the level where thepolarization changes. We denote this as a rotation from the x-ycoordinate frame, which is the initial coordinate frame of the sources,into x'-y' frame, the frame of the S-wave polarizations. The rotationsimulates lining up the x source polarization along the direction of thefast S-wave polarization of the upper layer. Ideally, after thisrotation, no signal energy would remain on the X'Y' or Y'X' componentsof the upper layer; and the signal on the Y'Y' components of the upperlayer should be time-lagged versions of the X'X' components.

After rotation into the primed coordinate frame comes the key step ofapplying a static shift to all data generated by one of the simulatedsource polarizations, the y', for example; thus the Y'X', Y'Y' and Y'Z'components from all depths at and below the bottom of the upper layerare time shifted by the amount needed to eliminate the lag between X'X'and Y'Y' wavelets at the bottom of the upper layer. Eliminating this lagis equivalent to positioning simulated x' and y' source polarizations atthe same depth, specifically at the top of the second layer. The initialrotation will not have properly minimized energy on the X'Y' or Y'X'components of the lower layer because the effective x' and y' sourcepolarizations acted as though they were excited at different depths(i.e., different times). The "rotations" which follow the stripping,however, should do a good job of minimizing energy on those off-diagonalcomponents down to the bottom of the second layer. Also, "rotations"after stripping should cause lags to increase from a value of zero atthe level where change occurs to progressively larger values. Of course,data will not ordinarily be recorded precisely where a change occurs, soeven in principle the lag should not always be strictly zero at thelevel closest to the interface.

The above described procedure for analyzing vertical seismic profileshear wave data, or surface seismic reflectors shear wave data may befurther described in the following manner. The data is defined to haveat least two linearly independent, nearly orthogonal, and nearlyhorizontal source axes. Each source axis has at least two correspondingreceiver axes.

1. An initial analysis of shear wave polarization directions relative toa fixed coordinate frame is performed, and apparent time lags betweenfast and slow shear waves are determined at several depths.

2. Cues in the data are identified that suggest shear wave polarizationchange.

3. The natural polarization directions of and the time lag between thefast and slow shear waves in an upper layer are determined, above andadjacent to the shallowest depth where the cues suggest polarizationchanges. Other depths may be used as well, even if no cues suggestpolarization changes.

4. The source and receiver axes of all the data that are below or at theshallowest depth of indicated polarization changes are then rotated byan azimuth angle determined down to this depth, so that the first sourceand receiver axes are aligned with the natural polarization direction ofthe fast shear wave, and the second receiver axis is at a significantlydifferent azimuth angle, and so that the second source and firstcorresponding receiver axis are aligned with the natural polarizationdirection of the slow shear wave in the upper layer, while the secondcorresponding receiver axis is at a significantly different azimuthangle.

5. A static shift is then applied to all data components correspondingto one of the effective sources, either to components corresponding tothe source aligned with the fast shear wave polarization direction or tocomponents corresponding to the source aligned with the slow shear wavepolarization directions, to eliminate the time lag in the upper layerabove and adjacent to the shallowest depth where the cues suggestpolarization changes are indicated.

The above method may also be used for vertical seismic profile (VSP)data or surface seismic reflection data that has only a single sourceaxis. Only the receiver axes are rotated in this case.

If surface seismic reflection shear wave data is analyzed, one methodincludes an initial analysis of shear wave polarization directionsrelative to a fixed coordinate frame in similarly recorded VSP data froma nearby well, and the subsequent determination of the time lags.Surface seismic reflection shear wave data can also be analyzed withoutthe use of VSP data.

EXAMPLE 1. Lost Hills

Data sets to be discussed in detail are from nine-component VSPsrecorded in two wells 862 ft apart, the 11-10X and the 1-9J wells of theLost Hills oil field in the southern San Joaquin valley of California.By nine-component data we mean records from three orthogonal receivercomponents which detected waves as if from three separate, orthogonalsource polarizations as illustrated on FIG. 3. The x-axis shown on FIG.3 is along a source vehicle axis, and receiver axes are computationallyrotated after recording to coincide with the source axes. The 2×2 S-wavedata matrix consists of four of the nine data components obtained withthree orthogonal sources and three orthogonal receivers. For example,the XY data component is from the x source component and the y receivercomponent. Except for preliminary processing, only data of the 2×2S-wave data matrix was treated; that is, data from x and y sources andreceivers, or four of the nine components. The coordinate frame forrecording and processing was a right-handed Cartesian frame with thex-axis along a source vehicle axis. After determining S-wavepolarization directions, we reoriented the frame relative to true north.

The 11-10X Well

For the 11-10X well, two orthogonally oriented Omnipulse airgun sourceswere used, and were located 57 and 68 ft from the well and as close toeach other as possible. Data were recorded without moving the sources.Source guns were tilted at 45°, and each was fired five times left andfive times right for a total of 20 pops per receiver level. Sourcezero-times were obtained from accelerometers screwed into thebaseplates. Locations and azimuths of sources were determined bysurveyors after we completed the VSP.

The downhole receiver was a three-component (3-C) SSC K tool with aGyrodata gyrocompass for determining absolute orientation. With thereceiver clamped at the maximum depth, 1720 ft, and sources at VSPpositions several series of source impacts were recorded before, duringand after the hydraulic fracturing of the 12-10 well to monitor anychanges in S-wave polarization that might result from the fracturing.Fracturing did not detectably affect data of the 11-10X well, althoughit caused transient changes in data simultaneously monitored in a wellopposite the 12-10 well. After recording at the 1720 ft depth, recordingoccured at increments of 80 ft from 1700 ft to 900 ft, with the finallevel at 800 ft.

The 1-9J Well

For the 1-9J well, a single ARIS (ARCO Impulsive Source, provided byWestern Geophysical) was used for the near offset portion of the VSP andalternated between two ARIS sources for rings of offset VSPs. For nearoffset recording, the ARIS was 50 ft from the well. For offsetrecording, we positioned the sources successively at eight pointsnominally 45° apart in each of two concentric rings nominally at 350 and700 ft from the well. Each source position was marked with two 14 inchrebar pegs whose locations were subsequently surveyed for accuratesource locations and azimuths. For near offset recording a special ARISbaseplate pad of riprap and road base gravel was built in order to doall recording without moving the source. For offset recording no padswere needed because source effort at a given position was small. ARISmade 20 impacts per receiver level or offset position, five in each offour directions--fore, aft, left and right--with the impactor tilted 15°from the vertical. The source vehicle axis pointed towards the well atevery source location. Source zero-times were obtained from pulses froman accelerometer atop the impactor, the pulses transmitted to therecording truck via hard-wire connection.

The downhole receiver for the 1-9J well was the LRS-1300 3-C tool withthe Gyrodata gyrocompass attached. Receiver components were gimballed sothat two were always horizontal and the third vertical. Recordingoccurred at increments of 100 ft from depths of 2100 ft to 100 ft forthe near offset VSP but at a fixed 2000 ft for the offset VSPs. Aftercompleting the near offset VSP the receiver was lowered to the 2000 ftlevel. That level was recorded again, without moving the source, beforegoing to the offset VSP locations.

Although the source baseplate was not moved for all near offsetrecording, it had sunk more than a foot from the beginning of nearoffset recording to the end. The receiver remained clamped at the 2000ft level without repositioning for all subsequent offset VSP recording.

The well was a nearly vertical cased and cemented hole which had not yetbeen perforated. Maximum deviation from vertical was 1.1°, and thebottom of the hole was laterally displaced only 10 ft from the top. Thefluid level was lowered to about 300 ft to avoid tube waves, which wereundetectable in both wells.

DATA CONDITIONING

Birefringence effects were analyzed in data which were thought to be asclose to being unprocessed as possible, but the following dataconditioning steps were deemed necessary. The first step was tocalculate and apply zero-time corrections (statics) based on sourceaccelerometer pulses. The second step was, for each receiver depth orsource offset position, and for each receiver component, to sum the fivetraces of like source polarity and then subtract sums for which impactswere azimuthally opposite in order to simulate a source that applied apurely horizontal impulse. Such a source produces vertically travelingS-waves with little contamination from vertically traveling P waves. Afurther conditioning step was to rotate the x-axis of the downholereceiver into alignment with the source axis, accomplished with the aidof gyrocompass and surveyor data. Also, data from receivers that werenot gimbal mounted were rotated initially to make the receiver z-axisvertical.

Before analyzing the data for S-wave polarization directions, wecomputationally rotated the receiver data so as to minimize S-waveenergy on the vertical components. This rotation, which requires twoEuler angles, causes the plane of the two nearly horizontal receivercomponents to coincide with the plane of S-wave displacements.

For near offset VSP data the amount of tilt was small, typically 6°-10°.Such a tilt puts the receiver plane out of alignment with the sourceplane, which is horizontal. However, this misalignment is unlikely tocause problems because of the small size of the tilt and because sourceradiation patterns put S-wave energy nearly equally into all possibleS-wave polarization directions for nearly vertical travel. Whether ornot the tilt was applied made no difference in azimuth angles and anegligible difference in lags (0.1 ms maximum) calculated from nearoffset data. For the offset data a few of the angles differed by 1°.

The final data conditioning steps involved amplitude adjustments andbandpass filtering. It is assumed that the components of body waves inthe y direction from the x oriented source must be identical to thecomponents of body waves in the x direction from the y oriented sourceThat is, to diagonalize the 2×2 S-wave data matrix by a single rotationangle, it is necessary that the XY and YX data components be identical,where XY indicates data from the x source on the y receiver. For thiscase of nearly vertical rays, and under the assumptions of nodifferential S-wave attenuation and isotropic geophone response, anydifferences in total wave energy from the x source relative to thosefrom the y source should be attributable to source or near surfaceproperties. Hence, an amplitude adjustment was applied to all data(i.e., data from all three receiver components) of the y source to makethem, in a time window corresponding to the S-wave wavelets, have thesame energy as those of the x source. For effectiveness of data displaythe energy of the data in that S-wave time window was also adjusted tobe the same at every depth, while taking care not to alter relativeamplitudes of data from a given source component. Finally, to eliminatehigh frequency noise, a mild high-cut filter was applied.

RESULTS Near Offset VSP Data

Data from the 1720 ft level of well 11-10X are shown in FIG. 4 asinitially recorded, and in FIG. 5 after "rotation" to minimize energy onthe off-diagonal components. The similarity of the two S-wave waveletsafter rotation is noteworthy, as are the relatively low amplitudes ofthe off-diagonal components.

Data from well 1-9J comprise a much more complete set than do those ofwell 11-10X. FIG. 6 shows the same data after rotation to minimizeenergy on the off-diagonal components in the analysis window indicated.Low amplitudes within the analysis window on the off-diagonal componentsat all depths suggest the rotation criterion worked well for this dataset, and that the subsurface S-wave polarizations were relativelyuniform.

Initial rotation analyses gave azimuths for the 11-10X well data thatwere nearly constant over their relatively limited depth range. Incontrast, analyses of the 1-9J data gave azimuths that showed asubstantial and systematic change with depth (FIG. 7.) The two datapoints in FIG. 7 at the 2000 ft level were recorded on successive days,one near the beginning of the experiment, the other after recording upto the 100 ft level and then lowering the tool again, and hence providea measure of reproducibility of results. The change in polarizationazimuth with depth was counter to expectations from models and led us tosuspect that a near surface layer with a different polarization azimuthwas contaminating analysis of deeper data. The strongest indication thatS-wave polarizations in the near surface layer were different from thoseat greater depths is the azimuth of 13° at 100 ft in FIG. 7. Subsequentangles show a systematic increase in azimuth angle, to 31° at 200 ft andfrom there up to about 60° at 2100 ft. The change is not erratic, asmight be expected for random errors, but smooth, indicating a possiblesystematic error that might be eliminated by stripping off a nearsurface layer.

Initial rotation analyses gave a time lag of 9.6 ms at 100 ft and 6.5 msat 200 ft, below which level the lag increased monotonically down to 800ft. The 9.6 ms value at 100 ft is assumed to be aberrant, possiblybecause of the horizontal component of raypath at that small depth.Consequently 6.5 ms was chosen as the amount to strip off initiallyafter rotating source and receiver x-axes to 13°, the azimuth calculatedat the 100 ft depth.

Stripping simplifies the picture considerably. Instead of S-wavepolarization aximuths, varying more than 25° for depths below 100 ft, asthey did before stripping, the azimuths now cluster tightly about 60°with standard deviation of 2.8°. The first five azimuths, however, showa systematic drop which look as suspicious as the previous systematicrise in azimuth. Hence it is suspected that the initial angle and lagare not optimal. To explore the dependence of the azimuths on initialangle and lag, the near surface layer is stripped off using severalother starting angles and lags. Results in FIGS. 8 and 9 show thatcalculated azimuth angles are insensitive to starting angle butsensitive to starting lag. An angle of 6° and a lag of 5 ms were chosenas the best values Comparing data analyses before and after strippingoff the near surface layer (FIG. 8) illustrates how a highlybirefringent, thin layer can contaminate analysis of data recorded morethan a 1000 ft below it.

Variations in S-wave lags with depth after stripping off the nearsurface layer (FIG. 10) indicate a significant change in birefringenceal about 700 ft. The lags rise uniformly, then level off, then dropbefore continuing to rise again. If the subsurface were homogeneous, thelags would continue to increase at a constant rate; while if the rockbecame isotropic, lags would remain constant. The only way lags candiminish, as they do from 900-1200 ft, is for anisotropy to change.

Stripping down to 700 ft and then performing the "rotation" analysisshowed that there was no significant change in azimuth and no consistentincrease in lags until below 900 ft. Azimuth changed between 900 ft and1200 ft, but changes in lags there were inconsistent and small, reachinga maximum of 2.1 ms. The zone from 900-1200 ft, then, caused too littleS-wave splitting to have a significant impact on polarization analysisbelow 1200 ft. The final layer stripping of the 1-9J well data henceinvolved stripping off the zone from 900-1200 ft. Note that the azimuthchange in this zone was undetectable before layer stripping (FIG. 7).

Results of layer stripping analyses from the surface to TD aresummarized in FIGS. 11 and 12. On FIG. 11, the angles and lags postedalongside the data points indicate values of layer stripping parametersapplied at the top of the layer. For example, the near surface layer wasstripped off with an initial rotation angle of either 7° or 0°,indicated by the different symbols, and a static of 5 ms. (These anglesunlike the others are relative to the source azimuth, which was N6°E.)Layer stripping parameters for deeper layers are given relative to theparameters for the layers immediately above them. The similarity of thetwo sets of results shows that a 7° difference in initial rotation anglehad little effect on answers at deeper levels. Except for the nearsurface layer and the zone from 700-1200 ft, the subsurface at the 1-9Jwell (FIG. 12) proved to be rather uniformly birefringent.

As a check on the validity of layer stripping, it is useful to monitorVSP traces closely after each stripping to determine whether the resultsfit layer stripping models. We have found that seismic data usually fitbetter after layer stripping than before. For example, FIG. 13 comparesoff-diagonal components at the deepest levels before and after layerstripping. Traces are from depths below 1200 ft. In the analysis window,signal amplitudes are lower after layer stripping (bottom traces) thanbefore (top traces). This indicates that layer stripping caused a betterfit with the seismic model. Wave amplitudes in the figure relative totrace spacing are four times those of FIG. 6. According to the model,amplitudes of the S-wave direct arrivals should be zero after the"rotations". Although FIG. 6 shows that amplitudes of off-diagonaldirect arrival S-waves are low relative to those of S-waves on thediagonals, they are clearly lower after layer stripping than before.

Offset VSP Data

Offset data of the 1-9J experiment gave remarkably consistent S-wavepolarization azimuths, the mean azimuth being 55° and the standarddeviation 6.3°. The consistency results from the high S/N, from therelative simplicity of anisotropy in that area and from the fact thatthe near surface layer had little influence on data recorded 2000 ftbelow it. The lags are much less consistent than the azimuths and varysystematically along the polarization direction of the fast S-wave. Itis likely that the variation in lags derives from shallow raypathsegments, because variations of the magnitudes indicated would beunlikely to originate from portions of raypaths in close proximity, suchas those at depth, which converge on the receiver.

Support for this proposed correlation between lag variations and shallowraypath segments comes from comparing lags of the 11-10X VSP with thoseof the near offset 1-9J VSP. The increase in 11-10X lags between1200-1700 ft resembles that of the 1-9J (FIG. 10). The absolutemagnitudes, however, are lower in the 1-9J data by about 8 ms,consistent with the lag variation observed. Part of the difference inabsolute lag (up to 4 ms) appears to result from a relatively smalldecrease in 11-10X lag in the anomalous zone from 900-1200 ft, but therest must occur shallower in the section.

LOST HILLS RESULTS

S-wave polarization azimuths are consistent with depth for a givenanisotropic layer at Lost Hills; that is, in the near offset VSP theyare consistent from 200-900 ft and from 1200-2100 ft. Consistency isnoteworthy because each calculated azimuth is the result of anindependent set of measurements. The lesser consistency in the deeperzone is expected, because layer stripping removes the "inertia" thatbuilds up in polarization determinations as the lag between the S-wavewavelets increases. The high overall consistency in polarization azimuthresults from consistency in subsurface properties, from highsignal/noise (S/N) in VSP direct arrivals and from the fact that wavesalong vertical raypaths satisfy the model assumptions employed in dataanalysis. The consistency in azimuth justifies the layer strippingmodel, which implicitly assumes that S-wave polarizations remainconstant over appreciable depth ranges, and suggests that we would haveobtained no better results by calculating transfer functions.

Layer stripping was effective and important for eliminating effects of athin, near surface anisotropic layer which had natural S-wavepolarizations different from those of deeper materials. Layer strippingwas less important for dealing with a change in anisotropy from 700-1200ft because of the small change in lag there. It is evident from dataanalysis (FIGS. 7 and 8) that the near surface layer adversely affectedpolarization analysis down at least to 1500 ft, and to a serious degreedown to about 600 ft; but the effect is small at the deepest levels.Birefringence of deeper formations will overcome contamination from anear surface layer when the lag from the near surface layer is smallcompared with a wavelength and when the lag between S₁ ' and S₂ ' ismuch larger than the lag from the near surface layer. FIG. 2 makes itapparent that, if the lag from the near surface layer were to approach asignificant fraction of a wavelength, layer stripping would be necessaryin order to analyze data by "rotations". However, when lags are thatlarge, other analysis techniques may work well.

While a preferred embodiment of the invention has been described andillustrated, it should be apparent that many modifications can be madethereto without departing from the spirit or scope of the invention.Accordingly, the invention is not limited by the foregoing description,but is only limited by the scope of the claims appended hereto.

What is claimed is:
 1. A method of analyzing seismic shear wave data,said data having at least one linearly independent, nearly horizontalsource axis (Sx), and each of said source axes having at least first andsecond linearly independent, nearly orthogonal, and nearly horizontalreceiver axes (Rx, Ry), to evaluate changes in shear wave polarizationwith depth comprising the steps of:(a) performing an initial analysis ofshear wave polarization directions relative to a fixed coordinate frame,and determining the apparent time lags between a fast and a slow shearwave, at several depths; (b) identifying cues in said data that indicateshear wave polarization changes; (c) determining the naturalpolarization directions of and the time lag between said fast and slowshear waves, in an upper layer, above and adjacent to the shallowestdepth where said cues suggest polarization changes are indicated; (d)rotating said receiver axes (Rx, Ry) of all of said data that are at orbelow the shallowest depth where said polarization changes areindicated, by an azimuth angle determined down to said depth, so thatsaid first receiver axis, (Rx) is aligned with the natural polarizationdirection of said fast shear wave in said upper layer, and said secondreceiver axis is at an azimuth angle which is really 90 degrees to saidfirst receiver axis; and (e) applying a static shift either to datacomponents aligned with said fast shear wave polarization direction, orto components corresponding to the source aligned with the slow shearwave polarization direction, to eliminate said time lag in said upperlayer above and adjacent to the shallowest depth where said cues suggestpolarization changes are indicated.
 2. The method of claim 1, furthercomprising the steps of:(f) repeating steps (a) and (b); and (g)repeating steps (c), (d), and (e) if step (b) identifies cues in saiddata that indicate further shear wave polarization changes.
 3. A methodof analyzing vertical seismic profile shear wave data, said data havingat least first and second linearly independent, nearly orthogonal, andnearly horizontal source axes (Sx, Sy), and each of said source axeshaving at least first and second linearly independent, nearlyorthogonal, and nearly horizontal corresponding receiver axes (Rxx, Rxy,Ryx, Ryy), to evaluate changes in shear wave polarization with depthcomprising the steps of:(a) performing an initial analysis of shear wavepolarization directions relative to a fixed coordinate frame, anddetermining the apparent time lags between a fast and a slow shear wave,at several depths; (b) identifying cues in said data that indicate shearwave polarization changes; (c) determining the natural polarizationdirections of and the time lag between said fast and slow split shearwaves, in an upper layer, above the adjacent to the shallowest depthwhere said cues suggest polarization changes are indicated; (d) rotatingsaid source and said receiver axes of all of said data that are at orbelow the shallowest depth where said polarization changes areindicated, by an azimuth angle determined down to said depth, so thatsaid first source axis (Sx) and said first corresponding receiver axis(Rxx) are aligned with the natural polarization direction of said fastshear wave in said upper layer, and said second corresponding receiveraxis (Rxy), is at an azimuth angle which is nearly 90 degrees to saidfirst receiver axis, and so that said second effective source axis (Sy)and said first corresponding receiver axis (Ryy) are aligned with thenatural polarization direction of said slow shear wave in said upperlayer, while said second corresponding receiver axis (Ryx) is at anazimuth angle which is nearly 90 degrees to said first receiver axis;and (e) applying a static shift to all data components corresponding toone of said effective sources, either to components corresponding to thesource aligned with said fast shear wave polarization direction, or tocomponents corresponding to the source aligned with the slow shear wavepolarization direction, to eliminate said time lag in said upper layerabove and adjacent to the shallowest depth where said cues suggestpolarization changes are indicated.
 4. The method of claim 3, furthercomprising the steps of:(f) repeating steps (a) and (b); and (g)repeating steps (c), (d), and (e) if step (b) identifies cues in saiddata that indicate further shear wave polarization changes.
 5. A methodof analyzing vertical seismic profile shear wave data, said data havinga single source axis, (Sx), oriented at an angle between naturalpolarization directions of said shear waves, and at least first andsecond linearly independent, nearly orthogonal, and nearly horizontalreceiver axis (Rx, Ry), to evaluate changes in shear wave polarizationwith depth comprising the steps of:(a) performing an initial analysis ofshear wave polarization directions relative to a fixed coordinate frame,and determining the apparent time lags between a fast and a slow shearwave, at several depths; (b) identifying cues in said data that indicateshear wave polarization changes; (c) determining the naturalpolarization directions of and the time lag between said fast and slowshear waves, in an upper layer, above and adjacent to the shallowestdepth where said cues suggest polarization changes are indicated; (d)rotating said receiver axes (Rx, Ry) of all of said data that are at orbelow the shallowest depth where said polarization changes areindicated, by an azimuth angle determined down to said depth, so thatsaid first receiver axis (Rx) is aligned with the natural polarizationdirection of said fast shear wave in said upper layer, and said secondreceiver axis (Ry), is at an azimuth angle which is nearly 90 degrees tosaid first receiver axis; and (e) applying a static shift either to datacomponents aligned with said fast shear wave polarization direction, orto components corresponding to the slow shear wave polarizationdirection, to eliminate said time lag in said upper layer above andadjacent to the shallowest depth where said cues suggest polarizationchanges are indicated.
 6. The method of claim 5, further comprising thesteps of:(f) repeating steps (a) and (b); and (g) repeating steps (c),(d), and (e) if step (b) identifies cues in said data that indicatefurther shear wave polarization changes.
 7. A method of analyzing andimproving surface seismic reflection shear wave data, said data havingat least first and second linearly independent, nearly orthogonal, andnearly horizontal source axes (Sx, Sy), and each of said source axeshaving at least first and second linearly independent, nearlyorthogonal, and nearly horizontal corresponding receiver axes (Rxx, Rxy,Ryx, Ryy), to evaluate changes in shear wave polarization with depthcomprising the steps of:(a) performing an initial analysis of shear wavepolarization directions relative to a fixed coordinate frame insimilarly recorded vertical seismic profile data from a nearby well, anddetermining the apparent time lags between a fast and a slow shear wave,at several depths; (b) identifying cues in said vertical seismic profiledata that indicate shear wave polarization changes; (c) determining thenatural polarization directions of and the time lag between said fastand slow shear waves, in an upper layer, above and adjacent to theshallowest depth where said cues suggest polarization changes areindicated in said vertical seismic profile data; (d) rotating saidsource and said receiver axes of all of said data that are at or belowthe shallowest depth where said polarization changes are indicated, byan azimuth angle determined down to said depth, so that said firstsource axis (Sx) and said first corresponding receiver ax is (Rx) arealigned with the natural polarization direction of said fast shear wavein said upper layer, and said second corresponding receiver axis (Rxy),is at an azimuth angle which is nearly 90 degrees to said first receiveraxis, and so that said second effective source axis (Sy) and said firstcorresponding receiver axis (Ryy) are aligned with the naturalpolarization direction of said slow shear wave in said upper layer,while said second corresponding receiver axis (Ryx) which corresponds tothat effective source an azimuth angle which is nearly 90 degrees tosaid first receiver axis; and (e) applying a static shift to all datacomponents corresponding to one of said effective sources, either tocomponents corresponding to the source aligned with said fast shear wavepolarization direction, or to components corresponding to the sourcealigned with the slow shear wave polarization direction, to eliminatesaid time lag in said upper layer above and adjacent to the shallowestdepth where said cues suggest polarization changes are indicted.
 8. Themethod of claim 7, further comprising the steps of:(f repeating steps(a) and (b); and (g) repeating steps (c), (d), and (e) if step (b)identifies cues in said data that indicate further shear wavepolarization changes.
 9. A method of analyzing surface seismicreflection shear wave data, said data having a single source axis, (Sx),oriented at an angle between natural polarization directions of saidshear waves, and at least first and second linearly independent, nearlyorthogonal, and nearly horizontal receiver axes (Rx, Ry), to evaluatechanges in shear wave polarization with depth comprising the stepsof:(a) performing an initial analysis of shear wave polarizationdirections relative to a fixed coordinate frame in similarly recordedvertical seismic profile data from a nearby well, and determining theapparent time lags between a fast and a slow shear wave, at severaldepths; (b) identifying cues in said vertical seismic profile data thatindicate shear wave polarization changes; (c) determining the naturalpolarization directions of and the time lag between said fast and slowsplit shear waves, in an upper layer, above and adjacent to theshallowest depth where said cues suggest polarization changes areindicated in said vertical seismic profile; (d) rotating said receiveraxes (Rx, Ry) of all of said data that are at or below the shallowestdepth where said polarization changes are indicated, by an azimuth angledetermined down to said depth, so that said first receiver axis (Rx) isaligned with the natural polarization direction of said fast shear wavein said upper layer, and said second receiver axis (Ry), is at anazimuth angle which is nearly 90 degrees to said first receiver axis;and (e) applying a static shift, either to data components aligned withsaid fast shear wave polarization direction, or to componentscorresponding to the slow shear wave polarization direction, toeliminate said time lag in said upper layer above and adjacent to theshallowest depth where said cues suggest polarization changes areindicated.
 10. The method of claim 9, further comprising the stepsof:(f) repeating steps (a) and (b); and (g) repeating steps (c), (d),and (e) if step (b) identifies cues in said data that indicate furthershear wave polarization changes.
 11. A method of analyzing and improvingsurface seismic reflection shear wave data, said data having at leastfirst and second linearly independent, nearly orthogonal, and nearlyhorizontal source axes (Sx, Sy), and each of said source axes having atleast first and second linearly independent, nearly orthogonal, andnearly horizontal corresponding receiver axes (Rxx, Rxy, Ryx, Ryy), toevaluate changes in shear wave polarization with depth comprising thesteps of:(a) performing an initial analysis of shear wave polarizationdirections relative to a fixed coordinate frame, and determining theapparent time lags between a fast and a slow shear wave, at severaldepths; (b) identifying cues in said data that indicate shear wavepolarization changes; (c) determining the natural polarizationdirections of and the time lag between said first and second shearwaves, in an upper layer, above and adjacent to the shallowest depthwhere said cues suggest polarization changes are indicated; (d) rotatingsaid source and said receiver axes of all of said data that are at orbelow the shallowest depth where said polarization changes areindicated, by an azimuth angle determined down to said depth, so thatsaid first source axis (Sx) and said first corresponding receiver axis(Rxx) are aligned with the natural polarization direction of said fastshear wave in said upper layer, and said second corresponding receiveraxis (Rxy), is at an azimuth angle which is nearly 90 degrees to saidfirst receiver axis, and so that said second effective source axis (Sy)and said first corresponding receiver axis (Ryy) are aligned with thenatural polarization direction of said slow shear wave in said upperlayer, while said second corresponding receiver axis (Ryx) is at anazimuth angle which is nearly 90 degrees to said first receiver axis;and (e) applying a static shift to all data components corresponding toone of said effective sources, either to components corresponding to thesource aligned with said fast shear wave polarization direction, or tocomponents corresponding to the source aligned with the slow shear wavepolarization direction, to eliminate said time lag in said upper layerabove and adjacent to the shallowest depth where said cues suggestpolarization changes are indicated.
 12. The method of claim 11, furthercomprising the steps of:(f) repeating steps (a) and (b); and (g)repeating steps (c), (d), and (e) if step (b) identifies cues in saiddata that indicate further shear wave polarization changes.
 13. A methodof analyzing surface seismic reflection shear wave data, said datahaving a single source axis, (Sx), oriented at an angle between naturalpolarization directions of said shear waves, and at least first andsecond inertly independent, nearly orthogonal, and nearly horizontalreceiver axis (Rx, Ry), to evaluate changes in shear wave polarizationwith depth comprising the steps of:(a) performing an initial analysis ofshear wave polarization directions relative to a fixed coordinate frame,and determining the apparent time lags between a fast and a slow shearwave, at several depths; (b) identifying cues in said data that indicateshear wave polarization changes; (c) determining the naturalpolarization directions of and the time lag between, said fast and slowshear waves, in an upper layer, above and adjacent to the shallowestdepth where said cues suggest polarization changes are indicated; (d)rotating said receiver axes (Rx, Ry) off all off said data that are ator below the shallowest depth where said polarization changes areindicated, by an azimuth angle determined down to said depth, so thatsaid first receiver axis (Rx) is aligned with the natural polarizationdirection of said fast shear wave in said upper layer, and said secondreceiver axis (Ry), is at an azimuth angle which is nearly 90 degrees tosaid first receiver axis; and (e) applying a static shift either to datacomponents aligned with said fast shear wave polarization direction, orto components corresponding to the slow shear wave polarizationdirection, to eliminate said time lag in said upper layer above andadjacent to the shallowest depth where said cues suggest polarizationchanges are indicated.
 14. The method of claim 13, further comprisingthe steps of:(f) repeating steps (a) and (b); and (g) repeating steps(c), (d), and (e) if step (b) identifies cues in said data that indicatefurther shear wave polarization changes.
 15. A method of analyzingseismic shear wave data, said data having at least one linearlyindependent, nearly horizontal source axis (Sx), and each of said sourceaxes having at least first and second linearly independent, nearlyorthogonal, and nearly horizontal receiver axes (Rx, Ry), to evaluatechanges in shear wave polarization with depth comprising the stepsof:(a) performing an initial analysis of shear wave polarizationdirections relative to a fixed coordinate frame, and determining theapparent time lags between a fast and a slow shear wave, at severaldepths; (b) determining the natural polarization directions of and thetime lag between said fast and slow shear waves, in an upper layer,above and adjacent to a selected depth where shear wave polarizationchanges are suspected; (c) rotating said receiver axes (Rx, Ry) of allof said data that are at or below the shallowest depth where saidpolarization changes are suspected, by an azimuth angle determined downto said depth, so that said first receiver axis, (Rx) is aligned withthe natural polarization direction of said fast shear wave in said upperlayer, and said second receiver axis is at an azimuth angle which isnearly 90 degrees to said first receiver axis; and (d) applying a staticshift either to data components aligned with said fast shear wavepolarization direction, or to components corresponding to the sourcealigned with the slow shear wave polarization direction, to eliminatesaid time lag in said upper layer above and adjacent to the shallowestdepth where said polarization changes are suspected.
 16. The method ofclaim 15, further comprising the steps of:(e) repeating step (a); and(f) repeating steps (b), (c), and (d) if further shear wave polarizationchanges are suspected.